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Summary 

A computer program has been written that allows the design of catalytic monoliths for closed-cycle 
carbon dioxide lasers. Using design parameters obtained from workers at NASA Langley 
Research Center and from the literature, several specific monoliths have been designed and the 
results have been communicated to the research group working on this project at Langley. Two 
oral presentations have been made at NASA-sponsored workshops - at Langley in January, 1988 
and in Gainesville, Florida in May, 1988. 

Progress Report 

A fully functional monolith design program has been written that inputs operating parameters and 
computes CO and O 2 conversions for specified a monolith length or, optionally, computes 
monolith length for specified reactant conversions. Currently, the program has more than 1,000 
lines of standard FORTRAN code and has been run on both a Macintosh personal computer and a 
VAX/VMS minicomputer. A program listing is given in the appendix of this report. 

Currently, a monolith geometry of parallel square channels is specified by the program. Variable 
thicknesses of active catalyst layers on the walls of the channels, as well as fully active walls, can 
be specified by the user. The program accounts for the following physical and chemical processes: 

• laminar flow of gas in the channels 

• mass transfer of reactants from the gas to the channel walls 

• diffusion of reactants in the porous catalyst layer 

• first-order-overall reaction kinetics in the porous catalyst layer 

• option for adiabatic or isothermal modes of operation of the monolith 

• heat release by the exothermic reaction in adiabatic operation mode 

• effect of temperature on the reaction rate 

• change in gas properties (viscosity, diffusivity, thermal conductivity) with temperature, 
pressure, and reactant conversion 

• calculation of gas pressure drop 
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As an example of the use of the design program, a design study will be presented here in which a 
specific monolith was designed. For a given set of conditions including monolith face dimensions, 
the optimal channel wall thickness and channel opening size were determined that minimize the 
monolith length and, thus, overall volume. 

Table 1 summarizes the parameters that were fixed and those that were varied in the design study. 
The face dimension of the monolith was fixed at 100 mm by 100 mm. For each design case, a 
wall thickness and a channel opening were chosen and the length of the monolith required to 
achieve 25% conversion of CO and O 2 at a pressure drop of 0.125 kPa (0.5 in. H 2 O) was 
computed. Figure 1 shows how the total number of channels and the channel opening changed as 
the wall thickness was varied. Each of the points in Figure 1 represent a specific monolith design 
that achieves the requirements of 25% conversion at a pressure drop of 0.125 kPa. As the wall 
thickness is increased and the number of channels in the monolith decrease, the size of the channel 
openings must increase in order to stay within the pressure drop constraint. 

Here, we use the term "% void" or void fraction to represent the volume of the total monolith that 
is taken up by the open channel passages. As shown in Figure 2, the void fraction of the monolith 
decreases as the wall thickness is increased. This decrease in void fraction means that the monolith 
tends to become more compact, that is, that the overall volume of the monolith decreases. The 
overall volume of the monolith, of course, is an important design criteria for aeronautical and space 
applications of closed-cycle carbon dioxide lasers. 

However, there is a factor that tends to offset this decrease in void fraction. As the wall thickness 
and the channel opening increase, resistances to mass transfer of reactants also increase: reactant 
molecules have greater distances to travel from the center of a monolith channel to the channel wall 
and from the channel wall into the porous catalyst layer. This effect increases the amount of active 
catalyst required to do the required conversion and is shown in Figure 3. 

Figure 4 shows how the two offsetting factors of decrease in void fraction and increase in mass 
transfer resistance couple to produce an optimal combination of wall thickness and channel opening 
that minimize the length and overall volume of the monolith. Figure 5 shows schematic 
representations of the relative dimensions of three of the design cases plotted in Figure 4: the left- 
most case (smallest wall thickness considered), the optimal case, and the right-most case in Figure 
4 (largest wall thickness considered). The optimal case monolith is only two-thirds of the length of 
the case considered with the thinnest walls, and is only slightly more heavy, as Figure 6 shows. 

Figure 7 compares the optimal design case considered here to a monolith with standard Coming 
Glass wall and channel dimensions. The optimal design is only about two-fifths the overall size of 
the monolith with off-the-self dimensions. 

Future Work 

The capabilities of the monolith design program will be expanded and the user interface will be 
enhanced. During the remainder of the project, the option to consider other channel geometries as 
well as parallel plate geometry will be added. In addition, the integration algorithm will be made 
more robust and efficient Further designs will be performed using NASA catalyst kinetics. 
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Table 1 

Design Study: 

Minimize Monolith Length, 

Keeping Pressure Drop Constant 


Square channels in square array 

Face dimension = 1 00 mm x 1 00 mm (4 in x 4 in) 

Gas flow rate = 1 liter/s (at inlet temperature of 300K) 

Porous wall material is uniformly active, with kinetics from 

Stark and Harris, J. Phys. E.: Sci. Instrum. 16, 492 (1983) 
Adiabatic operation 

CO and O 2 conversion = 25% (50 pmol O 2 /S for 0.5% O 2 in) 
Pressure drop across monolith = 0.125 kPa (0.5 in. H 2 O) 


Vary: 

Wall thickness 
Channel opening 

Minimize: 

Length of monolith 
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Total No. of 
Channels 


Figure 1 


Variation in No. Channels and Channel Opening 
at Constant Conversion and AP 
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% Void 


Figure 2 

Vo Void Decreases With Increasing Wall Thickness - 
Favoring More Compact Monolith 
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Figure 3 


Mass Transfer Resistance Increases with Wall Thickness and 

Channel Opening 
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Figure 4 

Minimize Length at Constant Conversion and AP 



1.38 mm wall, 

0.82 mm opening, 
14% void, 

209 mm long 
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Figure 5 

Minimize Length at Constant Conversion and AP 


0.548 mm opening 
0.27 mm wall 
45% void 
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305 mm 



10 mm 


0.82 mm opening 
1 .38 mm wall 
1 4% void 




209 mm 


10 mm 


1.18 mm opening 
3.0 mm wall 
8% void 



256 mm 
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Figure 6 


Slight Penalty in Monolith Weight Should be 
Offset by Savings in Housing Weight 



Weight (kg) 
at (0.55 g/cc) 
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Figure 7 

Optimum Design vs. Standard Corning Size 
(25 % Conversion) 


0.82 mm opening 
1 .38 mm wall 
1 4% void 



209 mm 

(AP = 0.5 in. H 2 0) 


10 mm 


1 .2 mm opening 
0.27 mm wall 
67% void 
(std. Coming size) 



502 mm 

(AP = 0.22 in. H 2 0) 
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APPENDIX: 

DRAFT OF COMPUTER PROGRAM LISTING, 6/30/88 
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C PROGRAM LASER17 WRITTEN BY KEITH GUINN AND SETH GOLDBLUM UNDER 
C THE DIRECTION OF PROF. RICHARD HERZ, CHEMICAL ENGINEERING, B-010 
C AT THE UNIVERSITY OF CALIFORNIA AT SAN DIEGO, LA JOLLA, CA 92093 
C FUNDED BY NASA LANGLEY RESEARCH CENTER, TECHNOLOGY 
C UTILIZATION AND APPLICATIONS OFFICE 
C MICROSOFT FORTRAN V2 . 2 MACINTOSH PLUS 
C JUNE 30, 1988 

C 

C PURPOSE: See introduction below 
C IDENTIFIERS: 

C A- RE ACT ION RATE CONSTANT OF CATALYST (mm A 3/g-cat s) 

C AMW- AVERAGE MOLECULAR WEIGHT OF GAS MIXTURE (g/mole) 

C AREF-VOL. PUMPING RATE OF CATALYST AT 2 98K (mm A 3/g-cat s) 

C ASK-QUERY DIRECTION NUMBER 

C AX (i, j) -i, j=l->6** CROSS TERM MATRIX FOR THERMAL CONDUCTIVITY 
C BEEP-BEEP CHARACTER 

C CGX-MONOLITH CHAN. 02 GAS CONC . AT POSITION X (mol/mm A 3) 

C CH - CHARACTER* 1 Dummy Variable 

C CHI- CHARACTER* 1 Dummy Variable 

C CHFLO-MONOLITH CHANNEL VOLUMETRIC FLOWRATE (mm A 3/S) 

C CHGVEL-MONOLITH CHANNEL GAS VELOCITY (mm/s) 

C CHHD-MONOLITH CHANNEL HYDRAULIC DIAMETER (mm) 

C CHSTV-MONOLITH CHANNEL SURFACE AREA TO VOLUME (1/mm) 

C CHWP -MONOLITH CHANNEL WETTED PERIMETER (mm,R) 

C CMCWT— CHARACTERISTIC MCWT(mm) 

C CONST1 , CONST2 , CONST3 , CONST4 , CONST5— RESULT HOLDING CONSTANTS 

C CONST5,CONST6,CONST7,CONST8, CONST 9-RESULT HOLDING CONSTANTS 

C CONSTA, CONSTB, CONSTC— RESULT HOLDING CONSTANT 

C CONV-CONVERS ION OF OXYGEN IN MONOLITH CHANNEL (%) 

C CONVX— CONVERSION OF OXYGEN AS FUNCTION OF X(%) 

C CP (i) , i=l, 6-HEAT CAPACITY OF SPECIES i AT TEMP. TGX (J/K-g) 

C CPGAS-HEAT CAPACITY OF BULK GAS MIXTURE (J/K-g) 

C CPCAT-HEAT CAPACITY OF CATALYST WASHCOAT (J/K-g) 

C CPSUP-HEAT CAPACITY OF MONOLITH SUPPORT (J/K-g) 

C CSubS-SUTHERLAND INTERACTION COEFFICIENT,! FOR NON-VERY POLAR 

C GASES 

C CSX— MONOLITH CHAN. WALL 02 GAS CONC AT POSIT. X (mol/mm A 3) 

C D (i) , i=l-6-BINARY DIFF. COEFF FOR SPECIES 02 & i,298K, 1 atm 
C (mm A 2/s) 

C DAB-DIFFUSIVITY OF 02 IN BULK GAS (mm A 2/S) 

C DABEFF-EFFECTIVE DIFFUSIVITY OF 02 IN WASHCOAT (mm A 2/S) 

C DABEFR-REF .EFFECTIVE DIFFUSIVITY OF 02 IN WASHCOAT, 2 98K 
C (mm A 2/S) 

C DCGX-CHANGE IN CGX FOR A CHANGE DX(mol/mm A 3) 

C DENGAS-DENSITY OF BULK GAS MIXTURE (g/mm A 3) 

C DENGSR- REFERENCE DENSITY OF BULK GAS MIXTURE (g/mm A 3) 

C DENCAT-DENSITY OF CATALYST WASHCOAT (g/mm A 3) 

C DENSUP-DENSITY OF MONOLITH SUPPORT (g/mm A 3) 

C DP-CHANGE IN PRESS ASSOC. WITH A CHANGE IN CHANNEL POSITION 
C DX (atm) 

C DTGX-CHANGE IN TGX FOR A CHANGE DX(mol/mm A 3) 

C DX- INCREMENTAL CHANGE IN MONOLITH CHANNEL POSITION X (mm) 

C DY5-CHANGE IN MOLE FRACT. OF 02 ASSOC. WITH A CHANGE DCGX IN 
C 02 CONC . 


NASA NAG-1-823, June 30, 1988 
(UCSD 88-5205) 


C ENGACT-CATALYST REACTION ACTIVATION ENERGY (J/mol) 

C FACTOR-PROFILE STEP SIZE FACTOR ( integer ) 

C FILNAM - NAME OF DATAFILE CONTAINING DEFAULT PARAMETERS 

C H-HEAT TRANSFER COEFFICIENT FROM GAS TO WALL ( J/mm^-S-K) 

C HEATRX- NEG. CHANGE IN ENTHALPY FOR 02+2CO->CO2 (J/mol) 

C I- INTEGER COUNTER 

C K-THERMAL GAS CONDUCTIVITY (J/mm-S-K) 

C J- INTEGER COUNTER 

C KC-MASS TRANSFER COEFFICIENT (mm/S) 

C KREF-REFERENCE THERMAL GAS CONDUCTIVITY ( J/mm^-S-K) 

C MCH-MONOLITH CHANNEL HEIGHT (mm) 

C MCST-MONOLITH CHANNEL SUPPORT WALL THICKNESS (mm) 

C MCW- MONOLITH CHANNEL WIDTH (mm) 

C MCWT-MONOLITH CHANNEL WASHCOAT THICKNESS (mm) 

C MFD-MONOLITH FACE DIMENSION (mm) 

C MON OL — COMP U TAT I ON LOOP COMPLETION MONOLITH LENGTH (mm) 

C MVFLO-MONOLITH VOLUMETRIC GAS FLOWRATE (mm A 3/S) 

C MW(i) — i->l, 6** MOLECULAR WEIGHT OF COMPONENT GASES (gms /mole) 

C N-EFFECTIVESS FACTOR FOR CATALYST WASHCOAL (NONDIM) 

C NOCHM-NUMBER OF CHANNELS IN MONOLITH (I) 

C NOCHV-NUMBER OF CHANNELS VERTICALLY IN MONOLITH 

C NOCHH-NUMBER OF CHANNELS HORIZONTALLY IN MONOLITH 
C 02 ZERO-MONOLITH INLET OXYGEN CONCENTRATION (%02) 

C P-DO LOOP VARIABLE 

C PHI (i, j) -i, j=l->6** CROSS TERM MATRIX FOR VISCOSITY 
C PRESS-PRESSURE OF GAS MIXTURE IN THE CHANNEL (atm) 

C PRESS 1-HOLDING CONSTANT FOR PRESSURE OF GAS MIXTURE 

C Q-HOLDING CONSTANT 

C R-GAS CONSTANT (J/mol-K) 

C Re - Reynold's number (dimensionless) 

C S-TAB CHARACTER 

C SS-CARRIAGE RETURM CHARACTER 

C SSS-A OR I FOR ADIABATIC OR ISOTHERMAL-CHARACTER 
C SSS1-P OR S FOR PROFILE OR SUMMARY-CHARACTER 

C SSS2-0 OR L FOR 02 CONVERSION OR MONOLITH LENGTH-CHARACTER 

C SU (i) -i->l, 6** SUTHERLAND CONSTANT OF COMPONENT 
C GASES (gms/mole) 

C SUINTER (i, j ) -i, j=l->6** SUTHERLAND INTERACTION MATRIX 
C T-GAS TEMPERATURE FOR VISCOSITY AND THERMAL CONDUCTIVITY 
C TB(i)-i->l, 6** BOILING POINT OF COMPONENT GASES (K) 

C TCON(i, j) -i->l, 6**, j->2, 7***, THERMAL COND . OF COMPONENT GASES 
C AT VARIOUS TEMPERATURES (mW/cm-K) 

C TCONT (i) -i->l, 6** THERMAL COND . OF COMPONENT GASES AT TEMP 
C T (mW/cm-K) 

C TCONDC-COMPOSITE GAS THERMAL CONDUCTIVITY (mW/cm-K) 

C TGX-MONOLITH CHANNEL BULK GAS TEMPERATURE AT POSIT. X (K) 

C THIELE-SLAB THIELE MODULUS FOR CATALYST WASHCOAT (NONDIM) 

C TREF-PHYSICAL PROPERTIES REF TEMP=298K (K) 

C TSX-MONOLITH CHANNEL WALL TEMPERATURE AT POSIT. X (K) 

C TZERO-MONOLITH INLET GAS TEMPERATURE (K) 

C VISC (i, j) -i->l, 6**, j->2, 7***, VISCOSITY OF COMPONENT GASES 
C AT VARIOUS TEMPERATURES (centipoises) 

C VISCT (i) -i->l, 6** VISCOSITY OF COMPONENT GASES AT TEMP 
C T (centipoises) 
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C VI SCTC— COMPOSITE GAS VISCOSITY (centipoises) 

C X-MONOLITH CHANNEL POSITION (mm) 

C Y (i) -i->l, 6 ** MOLE FRACTION OF COMPONENT GASES (gms/mole) 

C Y10,Y20,Y30,Y40,Y50, Y 60-HOLDING CONSTANTS FOR Y(i)'s 
C YT-SUM OF Y(i)-i->l ,6 

C **l=He, 2=Ar , 3=C02 , 4=CO, 5=02, 6=N2 
C ***2=200, 3=300, 4=400, 5=500, 6=600, 7=700K 
C ALGORITHM: 

C -Declare variables, open data file 

C -Introduction 

C -Initialize physical and option parameters to default values 
C -Display Selectable Parameter Summary 

C -Input physical parameters (if required) 

C -Monolith face dimensions 

C -Perform monolith parameter calculations 

C -Gas composition 

C -Inlet gas temperature, inlet gas flowrate 

C -Catalyst volumetric pumping rate and activation energy 

C -Choose operation options (if required) 

C -Adiabatic or isothermal operation 

C -Complete cone. & temperature profile output or summary 

C -Loop completion criteria (conversion or monolith length) 

C -Input desired final conversion or monolith length 

C -Computation loop step size, profile output step size. 

C -Allow revision of parameters if desired 

C -Compute viscosity and thermal conductivity parameters 

C -Read in molecular weight , diffusion coeff, thermal 
C conductivity, 

C viscosity, boiling point data 

C -Prepare to compute thermal conductivity & viscosity 

C -Print to screen and file, inputted, computed, and option 
C parameters 

C -Computation loop 

C -Intialize computation loop parameters 

C -Print profile header 

C -Step through computation loop until loop conditions are 

C met 

C -Print to data file and screen profile results at this 

C step 

C -Update physical properties that may change 

C -Query for further runs using different parameters 

C -Exit program 

q**********^ 0 q^j -0 vQ.;riclt) 3_0 s ************************** 

PROGRAM LASER17 

INTEGER NOCHV, NOCHH, NOCHM, FACTOR 

REAL* 8 MFD , MCW, MCH, MCWT, CMCWT, MCST, CHWP , CHSTV, CHHD 
REAL * 8 TZERO, 02ZER0, MVFLO, CONV, CHFLO, CHGVEL, MONOL 
REAL * 8 CGX, CSX, TGX, TSX, TREF, X, DX, DCGX, DTGX, AREF, A 
REAL* 8 KC, DENGAS, DENGSR, DENCAT, DENSUP, CPGAS, CPCAT, CPSUP 
REAL* 8 N, HEATRX, ENGACT, R, THIELE, DAB, DABR, DABEFF 
REAL* 8 H, K, KREF, CONST1 , CONST2 , CONST3 , P , Q, CONVX, ASK 

CHARACTER*! S, SS, SSS, SSS1 , SSS2 , BEEP, CH, CHI 
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CHARACTER* 20 FILNAM 

REAL * 8 T, Y(6) , YT, TCONDC, VISCTC, Y10, Y20, Y30, Y40, Y50, Y60 
REAL * 8 MW (6) ,TB(6) , SU(6) , CSub, DABEFR 
REAL *8 VISC (6, 2: 7) , SUINTER(6, 6) ,PHI (6, 6) 

REAL* 8 VISCT (6) ,TCONT(6) ,AX(6, 6) , AMW, PRESS, D ( 6 ) ,DY5,DP 

REAL * 8 CONST5, CONST4 , CONST6, CONST7, CONST8, CP ( 6) 

INTEGER I , J, CONSTA, CONSTB, CONSTC 

q* "k 'k * ~k ic introduct ion* ^'k^******************************'* - *********** 

Q***************************************************************** 

C Open Data file which will be used for output. 

OPEN ( UNIT=25, FILE= ' LDATA ' , STATUS= ’ NEW ' ) 

WRITE (9,*)' * * * * * *PROGRAM LASER17******* ' 

WRITE (9, *) 

WRITE (9,*)' The purpose of Laserl7 is to calulate the gas', 

2' concentration and temperature' 

WRITE (9,*) 'profiles for a monolith catalyst section of a', 

2' C02 laser. The C02 decomposes' 

WRITE (9, *)' when the laser is pulsed. The CO and 02', 

2' produced as a result of pulsing ' 

WRITE (9,*) 'are detrimental to the operation and efficiency', 

2' of the laser. The ' 

WRITE (9, *)' recovering rxn is CO+1/2 02 -> C02 . This program’, 
2' provides the means to model ' 

WRITE (9, *) 'the performance of a monolith catalyst section,', 
2' with various gas compositions,' 

WRITE (9, *) 'temperatures, catalyst activities, gas flowrates', 
2', oxygen conversion, monolith' 

WRITE (9, *)' face and length dimensions. Results can indicate', 
2' if constraints such as' 

WRITE ( 9, *)’ conversion, max. gas temperature, monolith weight', 
2' are satisfied and how ' 

WRITE (9,*) 'the system parameters may be altered to meet 
these ' , 

2' constraints. A number' 

WRITE (9,*) 'of parameters may be altered to meet these', 

2' constraints. A number' 

WRITE ( 9, *) ' of parameters can be supplied and a number of', 

2' options chosen to' 

WRITE ( 9, *)' customize the monolith design process. Default ', 

2 'values can also be used.' 

WRITE ( 9, * ) 'A review of the parameters and options chosen may', 
2 'be reviewed and altered' 

WRITE (9,*) 'prior to execution of the computation portion of', 
2' the program' 

WRITE (9, *) 

WRITE (9,*) ' (HIT RETURN TO CONTINUE) ’ 

PAUSE 
WRITE (9, *) 
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q "k -k ~k ★ ★ *Dg va. lues foir ★★★★ * ★★******★★★*★★★★★*★ ** 


First see if default parameters are already available in a data 
file. 

9000 WRITE (*,*) 

WRITE (*,*) 'Read in Data from Existing Data File? (Y/N) ' 

READ ( * , * ) CH 

IF ( (CH.EQ. ' Y' ) .OR. (CH.EQ. 'y ' ) ) THEN 

WRITE (*,*)' COOL .. . Thats Less Work for Both of Us.!!' 

WRITE (*, *) 

WRITE (*,*) 'What' s the name of the Existing Data file’ 

WRITE (*,*) ’to be opened?’ 

READ (*, *)FILNAM 

OPEN (UNIT=23, FILE=FILNAM, STATUS= ’ OLD ’ ) 

C If existing data file exists, read in data. 

READ (23, *) MFD 
READ (23, * ) MCST 
READ (23, *) MOW 
READ (23, *) MCH 
READ (23, * ) MCWT 
READ (23, *) Y (1) 

READ (23, * ) Y (2) 

READ (23, *) Y (3) 

READ ( 2 3 , * ) Y ( 4 ) 

READ ( 2 3 , * ) Y ( 5 ) 

READ (23, *) Y (6) 

READ (23, *) PRESS 
READ (23, *) SSS 
READ (23, *) SSS1 
READ (23, *) SSS2 
READ (23, * ) CONV 
READ (23, *) ENGACT 
READ (23, *) AREF 
READ (23, *) TZERO 
READ (23, * ) MVFLO 
READ (23, *) FACTOR 
READ (23, * ) DX 
CLOSE (UNIT = 23) 

C If no external data file exists, then assign default values. 
ELSE 

WRITE (*, *) 

WRITE (*,*)’ Default Values Being Assigned...' 

11 MFD=100 

MCST=0 . 5 
MCW=4 
MCH=4 
MCWT=0 .25 


1 6 


NASA NAG- 1-823, June 30, 1988 
(UCSD 88-5205) 


Y ( 1 ) =0 . 2 

Y (2) =0 . 0 
Y (3) =0 . 37 
Y ( 4 ) = 0 . 02 
Y{5)=0 .01 

Y ( 6 ) =0 . 40 
PRESS=1 . 0 
SSS= ' A ' 

SSS1='P ' 

SSS2= ' O 1 
CONV=5 

ENGACT=3 9700 
AREF=123 .4 
TZERO=300 
MVFLO=0 . 25 
FACTOR=5 
DX=1 . 0 

ENDIF 

C Compute a few needed parameters: 

Y10=Y ( 1 ) 

Y20=Y (2) 

Y30=Y (3) 

Y40=Y (4) 

Y50=Y (5) 

Y60=Y ( 6 ) 

PRESS 1=PRESS 
O2ZERO=100*Y(5) 

CHWP = 2 * ( MCW + MCH ) 

CHSTV = CHWP / ( MCW * MCH ) 

CHHD = 4 * ( MCW * MCH ) / CHWP 
NOCHV = INT (MFD / (MCW+2*MCST) ) 

NOCHH = NOCHV 
NOCHM = NOCHV * NOCHH 

0 ******** Print Out Main Manu **************************** 

9001 WRITE (*,*) 

WRITE (*, *) 

WRITE (*, *) 

WRITE (*, *) 

WRITE (*, *) 

WRITE (*, *) 

WRITE (*, *) 

WRITE (*, *) 

WRITE (*, *) 

WRITE (*, *) 

WRITE (*,*) 

WRITE (*, *) 

WRITE (*, *) 

WRITE (*,46) 
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46 FORMAT (1H ,T33,'Laser 17 Main Menu') 

WRITE (*, *) 

WRITE (*, *) 

WRITE (*,*) ' 1) Read In New Operating Parameters' 

WRITE (*,*) ' 2) Save Current Operating Parameters' 

WRITE (*,*) ' 3) Show Current Operating Parameters' 

WRITE (*,*) ' 4) Change Monolith Physical Dimensions' 

WRITE {*,*) ' 5) Change Monolith Inlet Parameters' 

WRITE {*,*) ' 6) Change Gas Flowrate and / or Temperature' 

WRITE (*,*) ' 7) Change Catalyst Reaction Rate Constant' 

WRITE (*,*) ' 8) Change Catalyst Activation Energy' 

WRITE ( * , * ) ' 9) Change Thermal Operation (Adia . /Isotherm) ' 

WRITE (*,*) ' 10) Change Output Profile (Full / Summary) ' 

WRITE (*, *) ' 11) Change Termination Parameter (Conv . /Len . ) ' 

WRITE (*,*) ' 12) Change Computation Loop Step Size’ 

WRITE (*,*) 

WRITE (*,*) ' 13) Compute Catalyst Profile' 

WRITE (*,*) ' 14) Exit Program' 

WRITE (*,*) 

WRITE ( *, * ) 'Type in Number Corresponding to Choice Above.' 

Read (*, *) ASK 

If (ASK.EQ.l) GO TO 9000 
If (ASK.EQ.2) GO TO 11000 
If (ASK.EQ.3) GO TO 9002 

If (ASK.EQ.4) GO TO 101 

If (ASK.EQ.5) GO TO 110 

If (ASK.EQ . 6) GO TO 114 

If (ASK.EQ. 7) GO TO 118 

If (ASK.EQ. 8) GO TO 119 

If (ASK.EQ. 9) GO TO 235 

If (ASK.EQ. 10) GO TO 236 

If (ASK.EQ. 11) GO TO 237 

If (ASK.EQ. 12) GO TO 309 

If (ASK.EQ. 13) GO TO 348 

If (ASK.EQ. 14) GO TO 99999 

WRITE (*,*) BEEP, BEEP, BEEP 
WRITE (*, *) 

WRITE (*,*) ' Invalid Choice, Ace.... ' 

WRITE (*,*) ' Try Again & Follow Directions.' 

WRITE (*,*) 

WRITE (*,*) 'Hit Return for Another Chance! ' 

PAUSE 
GOTO 9001 

q ★ * * * ★ ^O utput to screen for ^gvisw********************************** 
9002 Write (9, 6000) 

6000 FORMAT (20X, 'SELECTABLE PARAMETER SUMMARY’) 

WRITE (9, *) 

WRITE (9,*)' Monolith Dimensions (mm) :' 

WRITE (9, 6001) MCST*2 , MFD, MFD 

6001 FORMAT (' 1 . Support wall thickness : ' ,F6 . 2, 3X, 'Face dimension 


1 8 
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2F7 .2, ' x • ,F7 .2) 

WRITE (9, 6002) MCW, NOCHV, NOCHV 

6002 FORMAT (' 2 . Channel inner dimension :', F6 . 2, 3X, 'No .Face 
channels 1 , 

2 ' : ' , 14 , ' x ' , 14) 

WRITE (9, 6003) MCWT, 100* (MCW**2/ (MCST*2+MCW) **2) 

6003 FORMAT (' 3 . Washcoat thickness :',F6.2,' % open face',F6.1) 

WRITE (9,*) 

WRITE ( 9, * ) ' Monolith inlet parameters:' 

WRITE (9, 6004) Y(3),Y(4),Y(5),Y(6),Y(1),Y(2) 

6004 FORMAT (' 4 . Gas Composition (mole fraction): C02- ' , F5 . 3, 2X, 

2 ' CO -' ,F5.3,2X, '02 - ' , F5 . 3 , / , 35X, ' N2 - ' , F5 . 3, 2X, ' He - 

',F5.3,2X, 

3 ' Ar - ' , F5 . 3 ) 

WRITE (9, 6005) MVFLO, TZERO 

6005 FORMAT (' 5 . Gas Flowrate (1/s) F7 . 3, 9X, ' Gas Temperature (K) :' , 
2F7.2) 

WRITE (9, *) 

WRITE (9,*)' Catalyst Properties:' 

WRITE (9, 6006) AREF 

6006 FORMAT (' 6 . Reaction rate constant (mm A 3 02/gcat-s) : ' ,F8 .2) 

WRITE (9, 6007) ENGACT 

6007 FORMAT (' 7 . Activation energy (J/mol) :',F10.2) 

IF (SSS.EQ.'A') WRITE (9, 6008) 

6008 FORMAT (' 8 . Thermal Operation (adiabatic/isothermal) : Adiabatic ' ) 
IF (SSS.EQ.'I') WRITE (9, 6009) 

6009 FORMAT (' 8 . Thermal Operation 
(adiabatic/ isothermal) : Isothermal ' ) 

WRITE (9, *) 

WRITE (9, *) ' Computational loop parameters:' 

IF (SSS1.EQ. 'S') WRITE (9, 6010) 

6010 FORMAT (' 9 . Output file (Full Profile/Summary): Summary') 

IF (SSS1.EQ. 'P') WRITE ( 9, 6011) 

6011 FORMAT (' 9 . Output file (Full Profile/Summary): Full Profile') 

IF (SSS2.EQ. 'O') WRITE (9, 6012) CONV 

6012 FORMAT ('10. Termination on (02 conversion/Length):', 

2 '02 conversion. %02 : ' , F7.3) 

IF (SSS2.EQ. 'L') WRITE (9, 6013) MONOL 

6013 FORMAT ('10. Termination on (02 conversion/Length):', 

2 'Length. Length (mm) :' , F7.3) 

IF (SSS1.EQ. 'P') WRITE (9, 6015) DX, FACTOR 
6015 FORMAT ('ll. Computation loop step size (mm) = ',F6.3, 

2' Step size factor= ',14) 

IF (SSS1.EQ. 'S') WRITE (9, 60150) DX 
60150 FORMAT ('ll. Computation loop step size (mm) = ',F6.3) 

WRITE (*,*) 

WRITE (*,*) '==> Hit Return When Finished Viewing <==' 

PAUSE 
GOTO 9001 

QQ****************************************************************** 
Q * * * * * *P3.2T3.IT\Gt G Z* input SGCt ion************************************** 

101 WRITE (9, *) 
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WRITE (9,*) ' For a square face, input monolith face 
2 ' dimension (mm) ' 

READ (9,*) MFD 

IF (MFD. LE. 0.0) GO TO 101 

102 WRITE (9,*) 

WRITE (9,*)' Input the monolith channel support wall' 

WRITE (9,*)' thickness (mm) . This includes the washcoat . ' 

READ (9,*) MOST 
IF (MOST. LE. 0.0) GO TO 102 
MOST = MOST / 2.0 
IF (MCWT.GT.MCST) GO TO 104 

103 WRITE (9,*) 

WRITE(9,*)' For a square monolith channel, input monolith' 
WRITE (9,*) ' inner channel dimension (mm) . ' 

READ (9,*) MOW 

IF (MOW. LE. 0.0) GO TO 103 

MCH = MOW 

104 WRITE (9,*) 

WRITE (9,*)' Input monolith channel washcoat thickness (mm) . ' 
WRITE (9,*)' Washcoat cannot be >1/2 support wall thickness.' 
WRITE (9,*) ' Support wall thickness is', MCST*2,' (mm) ' 

READ (9,*) MCWT 

IF ( (MCWT. EQ.0) .OR. (MCWT.GT.MCST) ) GO TO 104 
C**********Perform monolith parameter calculations 
8900 CHWP = 2 * ( MOW + MCH ) 

CHSTV = CHWP / ( MCW * MCH ) 

CHHD = 4 * ( MCW * MCH ) / CHWP 

NOCHV = INT (MFD / (MCW+2 *MCST) ) 

NOCHH = NOCHV 
NOCHM = NOCHV * NOCHH 
GO TO 9001 

********Gas composition 

Define Mole fractions of species Y(i) 

The index of Y(i) correspond to the species. 

Species : l=He, 2=Ar , 3=C02, 4=CO, 5=02, 6=N2 . 

110 Write (9,*) 'Input section for gas composition' 

Write (9,*) 'The mixture can consist of only the following' 
Write (9,*) 'species: He, Ar, C02, CO, 02, N2 . The mole ' 

Write (9,*) 'fractions must add up to >=.99 to <=1.01. ' 
Write(9,*) 'All mole fractions must be > or = zero.' 

WRITE (9, *) 

Write (9,*) 'Input gas mole fraction for He' 

Read (9, *) Y(l) 

Write (9,*) 'Input gas mole fraction for Ar ' 

Read (9, *) Y(2) 

Write (9,*) 'Input gas mole fraction for C02 ' 

Read (9, *) Y(3) 

Write (9,*) 'Input gas mole fraction for CO' 

Read (9, *) Y(4) 

Write (9,*) 'Input gas mole fraction for 02' 

Read ( 9, *) Y(5) 

02ZER0=Y ( 5 ) *100 

Write (9,*) 'Input gas mole fraction for N2 ' 

Read (9,*) Y(6) 
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YT= Y ( 1 ) +Y (2 ) +Y (3) +Y(4) +Y(5) +Y(6) 

Y10=Y { 1) 

Y20=Y<2) 

Y30=Y (3) 

Y40=Y ( 4 ) 

Y50=Y (5) 

Y60=Y { 6) 

C Check if restrictions are met . If not, force user to reinput 
C values 

Do 112 1=1, 6 

IF ( (Y(I) .LT.O) .or. (Y(I) .GT.1.0) ) GOTO 110 
112 Continue 

IF ( (YT.GT.1.01) .or. (YT. LT.O. 99 )) GOTO 110 
GO TO 9001 

C********Inlet gas temperature 
114 WRITE (9,*) 

WRITE ( 9, *)' Input monolith inlet gas temperature (K) ' 

WRITE ( 9, *)' Temperature must be >200K and <700K.' 

READ (9,*) TZERO 

IF ( (TZERO. EQ. 0) .OR. (TZERO.LT. 200) ) GO TO 114 
IF (TZERO. GT. 700) GO TO 114 
C********Inlet gas flowrate 
116 WRITE (9,*) 

WRITE (9,*)' Input monolith inlet volumetric flowrate (1/S) ' 
READ (9,*) MVFLO 
IF (MVFLO. EQ.0) GO TO 116 
GO TO 9001 

C********Catalyst reaction rate constant and activation energy 

118 WRITE (9,*) 

WRITE (9, *) 'Catalyst reaction rate constant-2 98K (mm^3/gcat-s) ' 
READ (9,*) AREF 
IF (AREF.EQ.0) GOTO 118 
GO TO 9001 

119 WRITE(9,*) 

WRITE (9, *) 'Catalyst activation energy (J/mol) ' 

READ (9,*) ENGACT 
IF (ENGACT .EQ. 0) GOTO 119 
GO TO 9001 

C******Adiabatic or isothermal operation 

235 WRITE (9, *)' Isothermal (I) or Adiabatic (A) operation?' 

WRITE (9,*) 'Type the letter I or the letter A. HIT RETURN' 

READ (9,*) SSS 
IF (SSS.EQ. ' i ' ) SSS = 'I' 

IF (SSS.EQ. 'a' ) SSS = 'A' 

IF ( (SSS.NE. 'I') .AND. ( SSS.NE.’A')) GO TO 235 
GO TO 9001 

C ******* FULL CONCEN. AND TEMP. PROFILE (P) OR SUMMARY ONLY' 

236 WRITE ( 9, *) 

WRITE ( 9, *)' Full concentration and temperature profile (P) OR' 
2 ' summary (S) ? ' 

WRITE (9, *) 'Type the letter P or the letter S. HIT RETURN ’ 
READ (9,*) SSS1 
IF (SSS1.EQ. 'p') SSS1 = 'P' 

IF (SSS1.EQ. 's') SSS1 = 'S' 
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IF ( (SSS1.NE. 'P') .AND. ( SSS1 . NE . ’ S ' ) ) GO TO 236 
GO TO 9001 

C********COMPUT A T IO N LOOP COMPLETION CRITERIA 

237 WRITE (9,*) 

WRITE ( 9, *)' Computation loop completion criteria can be ' 
WRITE (9, *) ' chosen to depend on oxygen conversion (0) or' 
WRITE (9,*) 'monolith length (L) . ' 

WRITE (9, *) 'Type the letter 0 or the letter L. HIT RETURN.’ 

READ (9,*) SSS2 

IF (SSS2 .EQ. 'o' ) SSS2 = 'O' 

IF (SSS2 . EQ . ' 1 ' ) SSS2 = 'L' 

IF ( (SSS2 .NE. 'O' ) .AND . ( SSS2 . NE . ' L ' ) ) GO TO 237 
IF (SSS2 .EQ. 'O' ) THEN 

238 WRITE ( 9, *)' Input desired %oxygen conversion, 0 TO 100%' 
READ (9,*) CONV 

IF ( (CONV. GT. 100.0) .OR. (CONV. LT.0) ) GOTO 238 
END IF 

239 IF (SSS2 .EQ. 'L' ) THEN 

WRITE (9,*) 'Input desired monolith length (10 TO 500mm) ' 
READ (9,*) MONOL 
IF (MONOL. LT.0) GOTO 239 
IF (MONOL .GT. 500) GOTO 239 
ENDIF 

GO TO 9001 

c * * * * * *C0MP UT AT I0 N LOOP STEP SECTION 
309 WRITE (9,*) 

WRITE (9, *)' Input computation loop step size (mm).' 

WRITE (9, *) 'Maximum loop step size 1mm.' 

READ (9,*) DX 
IF (DX . GT .1.0) GO TO 309 
WRITE (9, *) 

320 IF (SSS1 .EQ. 'P ' ) THEN 
323 WRITE ( 9, * ) 

WRITE ( 9, *)' Input profile step size factor.' 

WRITE ( 9, *)' Prof ile step size factor must be >=1 . ' 

READ (9,*) CONST 1 
IF (CONST1 .LT. 1 .0) GOTO 323 
FACTOR=CONSTl/l 
ENDIF 
WRITE (9, *) 

GO TO 9001 

11000 WRITE (9, *) 

C Give user option to save present default parameters. 


WRITE (*, *) 

WRITE (*,*)' Would you like to save the data to a' 
WRITE (*,*)' NEW or EXISTING data file? (N/E) ' 
WRITE (*, *) 

READ (*, *) CHI 

IF ( (CHI .EQ. 'E' ) .OR. (CHI .EQ. 'e ’ ) ) THEN 
WRITE (*, *) 
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WRITE (*,*) 'What' s the name of the EXISTING' 

WRITE (*,*) 'Data file to be opened?’ 

READ (* ,*) FILNAM 

OPEN ( UNIT=24, F I LE=F I LNAM, STATUS= ' OLD ' ) 

ELSE 

WRITE (*, *) 

WRITE (*,*) 'What' s the name of the NEW' 

WRITE (*,*) 'Data file to be opened?' 

READ ( * , * ) FILNAM 

OPEN ( UNIT=24, FILE=FILNAM, STATUS= 'NEW' ) 

END IF 

C Write default parameters to data file. 

WRITE (24, * ) MFD 
WRITE (24, *) MOST 
WRITE (24, *) MOW 
WRITE (24, * ) MCH 
WRITE (24, *)MCWT 
WRITE (24, *) Y (1) 

WRITE (24, *) Y (2) 

WRITE (24, *) Y(3) 

WRITE (24, *) Y (4) 

WRITE (24, *) Y(5) 

WRITE (24, *) Y(6) 

WRITE (24, *) PRESS 
WRITE (24, *) SSS 
WRITE (24, *) SSS1 
WRITE (24, *) SSS 2 
WRITE (24, *) CONV 
WRITE (24, *) ENGACT 
WRITE (24, *) AREF 
WRITE (24, *) TZERO 
WRITE (24, *)MVFLO 
WRITE (24, *) FACTOR 
WRITE (24, *) DX 
CLOSE ( UNIT = 24) 

GOTO 9001 
348 CONTINUE 

**Computation of parameters needed for thermal cond. calcuation 
Read in molecular weights MW (I) into array. 

The index of MW correspond to the species . 

Species : l=He, 2=Ar, 3=C02 , 4=CO, 5=02 , 6=N2 . 

DATA MW ( 1 ) , MW (2 ) , MW (3 ) , MW ( 4 ) , MW (5) , MW (6) 

+ / 4.0026,39.948,44.01,28.0105,31.9988,28.0134/ 

C Read in binary diffusion coefficients D(I) into array. 

C The index of D correspond to the species in which 02 is 
C diffusiving 

C Note: There is no D(5) defined. These are reference values at 

C 298K, 

C 1 atm, units (mm / '2/s) 

C Species: l=He, 2=Ar, 3=C02, 4=CO, 6=N2 . 

DATA D(l) ,D (2) ,D(3) ,D(4) ,D(6) 
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+ /72. 9/21. 2/16. 4,15. 6,22. 5/ 

C Read in viscosity (VISC) values into matrix. 

C Units on VISC are centipoises . The values are read off a 
C nomagraph 

C in Perry’s Chemical Engineering Handbook, 5th ed., p3-210. 

C The first index of VISC corresponds to the species. 

C Species: l=He, 2=Ar, 3=C02 , 4=CO, 5=02 , 6=N2 . 

C The second index of VISC corresponds to the temperature. 

C Temperature: 2=200K, 3=300K, 4=400K, 5=500K, 6=600K, 7=700K 

DATA 

VISC (1,2) , VISC (1,3) , VISC (1,4) , VISC (1,5) , VISC (1,6) ,VISC(1,7) 

+/.014, .0188, .022, .027, .031, .034/ 

DATA 

VISC (2, 2) , VISC (2,3) ,VISC (2, 4) ,VISC(2,5) , VISC (2, 6) ,VISC(2,7) 
+/.0164, .0223, .0278, .0327, .0376, .041/ 

DATA 

VISC (3,2) , VISC (3, 3) ,VISC (3,4) , VISC (3,5) , VISC (3,6) , VISC (3,7) 
+/.0105, .0147, .019, .023, .027, .031/ 

DATA 

VISC (4, 2) ,VISC (4, 3) ,VISC(4,4) ,VISC(4,5) ,VISC(4,6) ,VISC(4,7) 
+/.0138, .018, .022, .0255, .029, .0325/ 

DATA 

VISC (5, 2) ,VISC (5,3) , VISC (5,4) ,VISC (5,5) , VISC (5, 6) , VISC (5, 7) 

+ /.0155, .02, .0247, .0288, .0327, .0368/ 

DATA 

VISC (6, 2) ,VISC (6, 3) ,VISC(6,4) ,VISC(6,5) , VISC (6, 6) ,VISC(6,7) 

+ /.013, .0175, .0218, .026, .029, .0332/ 

C Read in boiling point (TB) data into array. 

C The index of TB corresponds to the species. 

C Species: l=He, 2=Ar, 3=C02 , 4=CO, 5=02 , 6=N2 . 

C Units are K. The values are taken from "Properties of Gases 
C and Liquids", 3rd ed.,Reid, Prausnitz,& Sherwood, p629-633 
C and p509 (Note says to use 79K for He) . 

DATA TB ( 1 ) , TB ( 2 ) , TB ( 3 ) , TB ( 4 ) , TB ( 5 ) , TB ( 6 ) 

+/7 9, 87 .3, 194 .7, 145.2, 90 .2,77 .4/ 

C Compute Sutherland (SU) constant for each species 
C The index of SU corresponds to the species. 

C Species: l=He, 2=Ar, 3=C02, 4=CO, 5=02 , 6=N2 . 

C The Sutherland constant is defined on p509 of "Properties of 
C Gases 

C and Liquids", 3rd ed.,Reid, Prausnitz,& Sherwood. 

SU (i) =1 . 5*TB (i) . 

DO 150 1=1, 6 

SU ( I) =1 . 5*TB (I) 

150 CONTINUE 

C Compute the Sutherland interaction constants, SUINTER 

C The indexes of SUINTER corresponds to the species. 

C Species: l=He, 2=Ar, 3=C02, 4=CO, 5=02, 6=N2 . 

C The Sutherland interaction constant is defined on p509 of 

C "Properties of Gases and Liquids", 3rd ed.,Reid, Prausnitz, 

C & Sherwood. SU (i) =1 . 5TB (i) . None of the species are "very 
C polar, 
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C therefore, CSubS=1.0, see p509-510. 

CSubS=l .0 
DO 175 1=1, 6 

DO 170 J=l, 6 

SUINTER(I, J)=CSubS* (SU(I) *SU(J) ) **0.5 
170 CONTINUE 

175 CONTINUE 

C COMPUTE ADDITIONAL MONOLITH PARAMETERS 

C CONVERT MVFLO FROM 1/S TO mm A 3/S 

MVFLO = MVFLO *1.0E+06 
CHFLO = MVFLO / NOCHM 
CHGVEL = CHFLO / ( MCW * MCH) 

C**Output of selectable parameter to data file********************* 

WRITE (9,*)' *****SEE FILE LDATA FOR RESULTS*****' 

WRITE (25, 6100) 

6100 FORMAT (2 OX, 'SELECTABLE PARAMETER SUMMARY') 

WRITE (25,*)' Monolith Dimensions (mm) :' 

WRITE (25, 6101) MCST*2, MFD, MFD 

6101 FORMAT (' 1 . Support wall thickness : ' ,F6 . 2, 3X, 'Face dimension 

. t 

2F7.2, ' x ' , F7 . 2) 

WRITE (25, 6102) MCW, NOCHV, NOCHV 

6102 FORMAT ( ' 2 . Channel inner dimension :' ,F6 . 2, 3X, 'No. Face 
channels ' , 

2' : ',14, ' x ’,14) 

WRITE (25, 6103) MCWT, 100* (MCW**2/ (MCST*2+MCW) **2) 

6103 FORMAT ( ' 3 . Washcoat thickness :’,F6.2,' % open face’,F6.1) 

WRITE (25, *) 

WRITE (25,*)' Monolith inlet parameters:' 

WRITE (25, 6104) Y(3),Y(4),Y(5),Y(6),Y(1),Y(2) 

6104 FORMAT (' 4 . Gas Composition (mole fraction): C02- ' , F5 . 3, 2X, 

2 ' CO - ' , F5 . 3 , 2X, ' 02 - ' , F5 . 3, /, 35X, ' N2 - ' , F5 . 3, 2X, ’ He - 

' , F5 . 3, 2X, 

3'Ar - ' , F5 . 3 ) 

WRITE (25, 6105) MVFLO/1 . 0E+06, TZERO 

6105 FORMAT (' 5 . Gas Flowrate (1/s) F5 . 3, 9X, ' Gas Temperature (K) :' , 
2F7.2) 

WRITE (25, *) 

WRITE (25,*)' Catalyst Properties:' 

WRITE (25, 6106) AREF 

6106 FORMAT (' 6 . Reaction rate constant (mm A 3 02/gcat-s) : ' ,F8 .2) 

WRITE (25, 6107) ENGACT 

6107 FORMAT (' 7 . Activation energy (J/mol) :',F10.2) 

IF (SSS.EQ.'A') WRITE (25, 6108) 

6108 FORMAT (' 8 . Thermal Operation (adiabat ic/isothermal) : Adiabatic ' ) 
IF (SSS.EQ.'I') WRITE (25, 6109) 

6109 FORMAT (' 8 . Thermal Operation 
(adiabatic/isothermal) : Isothermal ' ) 

WRITE (25, *) 

WRITE(25,*)' Computational loop parameters:' 

IF (SSS1.EQ. ’S') WRITE(25, 6110) 

6110 FORMAT (' 9 . Output file (Full Profile/Summary): Summary') 
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IF (SSS1.EQ. 'P') WRITE (25, 6111) 

6111 FORMAT (' 9 . Output file (Full Profile/Summary): Full Profile') 

IF (SSS2 .EQ. 'O' ) WRITE (25, 6112) CONV 

6112 FORMAT ('10. Termination on (02 conversion/Length):', 

2 '02 conversion. %02 : ' , F7.3) 

IF (SSS2.EQ. 'L') WRITE (25, 6113) MONOL 

6113 FORMAT ('10. Termination on (02 conversion/Length):', 

2 'Length. Length (mm) :' , F7.3) 

WRITE (25, 6115) DX, FACTOR 

6115 FORMAT ('ll. Computation loop step size (mm) = ',F6.3, 

2' Step size factor= ’,14) 

6114 WRITE (25,*) 

C****** IN TIALIZE PHYSICAL PROPERTIES PARAMETERS******************** 

Q***ic****j<'k-k'k'k**'kic-k-k-k-k-k1c-k-k-k’k-k1<:-k-k'k'k'k'k'k'k'k'k'k'k’k’k'k'k'k'k’k*'k'k'k'k , k , k‘k’k'k’k'kic'k'k'k'k’k 

HEATRX=5 . 64E+05 
TREF=2 98 
DENCAT=5 . 0E— 04 
DABEFR=1 . 0 
R=8 .3144 

C*******EULER'S METHOD COMPUTATION OF CG (X) , TG (X) , CS (X) , &TS (X) ***** 

c****************************************************************** 

TGX= TZERO 

CGX= 4.092E-08 * 02ZERO /100.0 * (298/TGX) 

CONST3=CGX 

CGXEND =CGX * ( 1 - CONV/100) 

S = CHAR (9) 

SS = CHAR ( 13 ) 

BEEP= CHAR (7) 

X=0 

CMCWT= ( (MCH+2 *MCWT) * (MCW+2*MCWT) - (MCH*MCW) ) /2/ (MCH+MCW) 

WRITE (9, 9950) 

WRITE (25, 9950) 

9950 FORMAT ( ' X-mm ' , 2X, ' %CONVER ’ , 2X, ' 02gas<-mM0L/L->02surf ', 
+' Tgas<-K->Tsurf ETA DP (KPA) ') 

c ********** LOOP UNTIL REQUIRED CONDITION IS MET************ 
0 ************************************************************ 
CONSTA=0 

DO 8000 P=0, 1000, DX 
IF (SSS2 .EQ. 'O' ) THEN 

IF (CGX .LT. (CGXEND / (TZERO/TGX) ) ) GO TO 9900 
GO TO 8250 
ENDIF 

IF (SSS2.GE. 'L') THEN 

IF (P.GT. MONOL) GO TO 16900 
IF (CGX. LE. 0.0) GO TO 16900 
ENDIF 

8250 AMW=0 

DO 8350 1=1, 6 

AMW=MW ( I ) * Y ( I ) +AMW 
8350 CONTINUE 

DENGAS=AMW* (PRESS/9896) /R/TGX 
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DAB=0 

DO 8450 1=1,6 

IF (I.EQ.5) GO TO 8450 
DAB=Y (I) / (D (I) * (TGX/298) **1 . 5) +DAB 
8450 CONTINUE 

DAB= (1-Y (5) ) /DAB/PRESS 
DABEFF=DABEFR* (TGX/298) **1 . 5 

C ********* COMPUTE THERMAL CONDUCTIVITY********************** 

IF (MOD ( CONSTA, 2) .GT. 0.0 ) GO TO 6666 
T=TGX 

C Compute viscosity (VISCT(i)) for a given species at temp. T 
C Units on VISCT are centipoises. 

C The index of VISC corresponds to the species. 

C Species: l=He, 2=Ar, 3=C02, 4=CO, 5=02 , 6=N2 . 

C CONST5 is a temporary value holding constant 
DO 200 1=1,6 

CONST5=VISC (I, INT (l+T/100 . 0) ) -VISC (I, INT (T/100 .0) ) 

CONST 5=CONST5* (T/100 .0- INT (T/100) ) 

VISCT (I) =VISC (I, INT (T/100 . 0) ) . + CONST5 
200 CONTINUE 

C Compute thermal cond. (TCONT (i) ) for a given species at temp. T 
C Units on TCON are mW/cM-K. 

C The index of TCONT (i) corresponds to the species. 

C Species: l=He, 2=Ar, 3=C02, 4=CO, 5=02 , 6=N2 . 

C The following polynomial equations used to compute the 
C thermal conductivity of each species at temperature T 

C have been taken from page 515 of "Properties of Gases and 

C Liquids" 

C 4th Ed. Reid, Prausnitz, & Sherwood 

C CONST5 is a temporary value holding constant 

TCONT ( 1) = (3.722E-2 + 3.896E-4*T - 7 . 450E-8*T**2 + 

1 1.290E-11*T**3)*10.0 

TCONT ( 2 ) = (2.714E-3 + 5.540E-5*T - 2 . 178E-8*T**2 + 

1 5 . 528E-12*T**3) *10.0 

TCONT ( 3 ) = (-7.215E-3 + 8.015E-5*T + 5 . 477E-9*T**2 - 
1 1 . 053E-11*T**3) *10 .0 

TCONT ( 4 ) = (5.067E-4 + 9.125E-5*T - 3 . 524E-8*T**2 + 

1 8 . 199E-12*T**3) *10.0 

TCONT (5) = (-3.273E-4 + 9.966E-5*T - 3 . 743E-8*T**2 + 

1 9.732E-12*T**3) *10.0 

TCONT ( 6 ) = (3.919E-4 + 9.816E-5*T - 5 . 067E-8*T**2 + 

1 1 . 504E-11*T**3) *10.0 


C Compute parameter matrix AX(i,j) 

C The indexes of A(i,j) correspond to the species. 

C Species: l=He, 2=Ar, 3=C02, 4=CO, 5=02, 6=N2 . 

C The A(i,j) matrix is defined on p509 of "Properties of Gases 
C and 
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C Liquids", 3rd ed.,Reid, Prausnitz, & Sherwood. 

C CONST5 is a temporary value holding constant 

CONST5=0 
DO 260 1=1, 6 

DO 270 J=l, 6 

CONST 5= (T+SU (I) ) / (T+SU ( J) ) * (MW(J) /MW (I) ) **0 .75 
CONST5= (1+ (CONST5 *VISCT ( I ) /VISCT(J) ) **0.5) **2.0 
AX (I, J) =CONST5/4 . 0* (T+SUINTER (I, J) )/(T+SU(I) ) 
270 CONTINUE 

260 CONTINUE 

C Compute the thermal conductivity of the mixture at TGX 
C CONST4 , CONST5 are holding constants 
CONST4=0 
DO 550 1=1, 6 
CONST5=0 
DO 650 J=l, 6 

CONST5=CONST5 + Y(J)*AX(I,J) 

' 650 CONTINUE 

CONST4=CONST4 + Y ( I) *TCONT ( I ) /CONST5 
550 CONTINUE 

TCONDC=CONST4 
K=TCONDC/ 10000.0 

C Compute parameter matrix PHI(i,j) 

C The indexes of PHI(i,j) correspond to the species. 

C Species: l=He, 2=Ar, 3=C02 , 4=CO, 5=02 , 6=N2 . 

C The PHI(i,j) matrix for Wilke's Approximation method 
C is defined on p411 of "Properties of Gases and 
C Liquids", 3rd ed.,Reid, Prausnitz, & Sherwood. 

C CONST5 is a temporary value holding constant 

DO 280 1=1, 6 
DO 290 J=l, 6 

CONST5= (VISCT ( I) /VISCT(J) ) **0.5* (MW(J) /MW(I) ) **0.25 
PHI (I, J)= (1+CONST5) **2.0/ (8* (1+MW(I) /MW ( J) ) ) **0.5 
290 CONTINUE 
280 CONTINUE 

C Compute the viscosity of the mixture at TGX 
C CONST4 , CONST5 are holding constants 
CONST4=0 
DO 750 1=1, 6 
CONST5=0 
DO 850 J=l, 6 

CONST5=CONST5 + Y ( J) *PHI (I, J) 

850 CONTINUE 

CONST4=CONST4 + Y (I) *VISCT (I) /CONST5 
750 CONTINUE 

VISCTC=CONST4 

6666 KC = 3.66 * DAB / CHHD 
H = 3.66 * K / CHHD 

C ********* C0 MPUTE heat capacity of gas mixture AT TGX******** 
£*************★**********★**★**★************★**************** 
CP (1) =20 . 786 
CP (2) =20 .786 

CP (3) =20.712+6. 2501E-02 *TGX-1 . 188E-5*TGX*TGX 
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+ -5 . 5773E-08*TGX*TGX*TGX + 4 . 1875E-11*TGX*TGX*TGX*TGX 

CP (4) =28. 777 + 5 . 5297E-03*TGX -3 . 1851E-5*TGX*TGX 

+ +7 . 094 6E-08*TGX*TGX*TGX -4 . 2501E-11*TGX*TGX*TGX*TGX 

CP (5) =31 . 033 - 2 . 0261E-02*TGX +6 . 3410E-5*TGX*TGX 

+ -4 . 9951E-08*TGX*TGX*TGX +9 . 1266E-12*TGX*TGX*TGX*TGX 

CP (6) =28 . 567 + 7.397 8E-03*TGX -3 . 5831E-05*TGX*TGX 

+ +6. 9415E-08*TGX*TGX*TGX -3 . 8499E-11*TGX*TGX*TGX*TGX 

CPGAS = Y ( 1 ) *CP ( 1 ) /MW ( 1 ) + Y (2) *CP (2) /MW (2) + 

+ Y (3) *CP (3) /MW (3) + Y (4) *CP (4) /MW(4) + 

+ Y (5) *CP (5) /MW (5) + Y (6) *CP (6) /MW(6) 

A = AREF * EXP ( -ENGACT / R * ( 1 / TGX -1 / TREF) ) 

THIELE = CMCWT * ( A * DENCAT / DABEFF) ** 0.5 

N = TANH ( THIELE ) / THIELE 

C0NST1= N*A*CMCWT*DENCAT 

CSX = CGX * { KC) / (KC+CONST1) 

TSX =TGX + HE ATRX*CSX* CONST 1/H 

DCGX = DX*KC*CHSTV*MCH*MCW/ CHFLO* (CSX-CGX) 

IF (SSS.EQ.’I') DTGX=0 
IF (SSS.EQ.'A') THEN 

DTGX = DX*H* CHSTV*MCH*MCW/ CHFLO* (TSX-TGX) /DENGAS/ CPGAS 
END IF 

CONVX = 100 * ( CONST3-CGX* (TZERO/TGX) ) /CONST3 
WRITE RESULTS TO SCREEN 

IF (MOD ( CONSTA, FACTOR) . GT . 0.0) GO TO 7777 

WRITE (25, 10) X, S , CONVX, S, CGX*lE+0 9, S, CSX* IE +0 9, S, TGX, S,TSX, S,N 
+, S, (PRESS-PRESS1) *101 .33 

FORMAT ( ' ' ,F6. 1, Al, 1X,F5.2, Al, 2 (3X,F7 ,4,A1) ,2 (3X,F8 .3, Al) ,F6.4 

+,A1,F7 .3) 

WRITE (9, 10) X,S, CONVX, S,CGX*1E+09,S,CSX*1E+09,S, TGX, S, TSX, S,N 
+,S, (PRESS-PRESS1) *101 .33 

CGX = (CGX +DCGX) * (TGX/ (TGX+DTGX) ) 

TGX = TGX +DTGX 
X = X + DX 

DP = 32* ( VISCTC/1000 ) *CHGVEL*DX/CHHD/CHHD/101330 . 0 
PRESS = PRESS - DP 

CONST7 = Y (1) +Y (2) +Y (3) +Y (4) +Y (5) +Y (6) 

DY5 = DCGX*R*TGX/ (PRESS/9869) 

CONST8 = Y ( 1 ) +Y ( 2 ) +Y (3 ) +Y (4) +Y (5) +Y (6) -3*DY5 

Y ( 1 ) = Y(l) *CONST7/CONST8 

Y (2) = Y (2) *CONST7/CONST8 
Y (3) = Y(3) *CONST7/CONST8 

Y ( 4 ) = (Y (4) -2*DY5) *CONST7/CONST8 

Y (5) = (Y (5) -DY5) *CONST7/CONST8 

Y ( 6) = Y(6) *CONST7/CONST8 
CONSTA=CONSTA + 1 

Check to make sure that Temperature, Pressure, and Reynold's 
number are still in valid range. 

Temperature must be between 200-700 K 
Pressure must be between 10E-03 and 3.45 bar 
Reynolds number must be in laminar region. < 2300 

RE=1000 . 0*CHHD*CHGVEL*DENGAS/VISCTC 
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IF ( (RE. GT. 2300.0) .OR. (TGX . GT . 700 . ) .OR. (TGX . LT . 200 . ) 

1 .OR. (PRESS. LT. 0.0.01) .OR. (PRESS. GT. 3.45) ) THEN 

WRITE (*, *) 

WRITE (*,*) 'Error : Parameters have exceeded' 

WRITE (*,*)' allowable Range... ' 

WRITE (*,*) 1 Returning to Main Menu' 

WRITE (*, *) 

WRITE (*,*) 'Hit Return When Finished Viewing' 

PAUSE 
GOTO 9001 

END IF 

8000 CONTINUE 

9900 WRITE (25, 10) X, S, CONVX, S, CGX*lE+0 9, S, CSX*lE + 0 9, S, TGX, S, TSX, S, N 
+ , S, (PRESS-PRESS1) *101.33 

WRITE (9,10) X,S, CONVX, S,CGX*1E+09,S,CSX*1E + 0 9,S, TGX, S, TSX, S,N 
+, S, (PRESS-PRESS1) *101 . 33 
16900 WRITE (25, *) 

WRITE (9,*) BEEP, BEEP, BEEP 
17000 WRITE (*,*) 

WRITE (*,*) ' ==> Hit Return to Return to Main Menu <==' 

PAUSE 

MVFLO=MVFLO/lE+ 6 
Y (1)=Y10 
Y (2 ) =Y20 
Y (3 ) =Y30 
Y (4 ) =Y40 
Y(5)=Y50 
Y ( 6) =*Y60 
PRESS=PRESS1 
GOTO 9001 

99999 END 
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